Transport properties for liquid silicon-oxygen-iron mixtures at Earth's core conditions 
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We report on the thermal and electrical conductivities ol two liquid silicon-oxygen-iron mixtures 
(Feo.82Sio.10Oo. 08 and Feo.79Sio.08Oo. 13), representative of the composition of the Earth's outer core 
at the relevant pressure-temperature conditions, obtained from density functional theory calculations 
with the Kubo-Greenwood formulation. We find thermal conductivities k =100 (160) W m~^ K~^, 
and electrical conductivities a = 1.1(1.3) x 1Q^Q,~^ at the top (bottom) of the outer core. 

These new values are between 2 and 3 times higher than previous estimates, and have profound 
implications for our understanding of the Earth's thermal history and the functioning of the Earth's 
magnetic field, including rapid cooling rate for the whole core or high level of radiogenic elements in 
the core. We also show results for a number of structural and dynamic properties of the mixtures, 
including the partial radial distribution functions, mean square displacements, viscosities and speeds 
of sound. 

PACS numbers: 



I. INTRODUCTION 

Transport properties of the Earth's core are of great importance to understand the thermal and magnetic behaviour 
of our planet. It is widely believed that the Earth's core is a mixture of iron and light impurities [T]. First principles 
calculations indicate that the liquid outer core contains a significant fraction of oxygen (between 8 or 13 %), and silicon 
and/or sulphur (between 8 and 10 %) [2]. Oxygen is expelled into the liquid upon freezing of the inner core at inner- 
core boundary (ICB), which helps to drive convection responsible for the generation of the Earth's magnetic field. 
Vigorous convection keeps the outer core in a well-mixed state, with a temperature distribution which closely follows 
an adiabatic profile. Heat flows from the bottom to the top of the core by thermal conduction and by convection. 
High thermal conductivity reduces the heat available to drive convection (in the limit of infinite thermal conductivity 
the whole core would become isothermal and would not convert). 

The electrical resistivity - and hence its inverse, the electrical conductivity a - determines the magnitude of Ohmic 
losses of the electric currents in the outer core that generate the magnetic field, and also the magnetic diffusion time 
(that is, the time that it would take for the magnetic field to decay in the absence of a generating mechanism). First 
experimental measurements of the electrical conductivity of pure iron combining high pressures and temperatures 
date back to the late sixties. They were performed on shock-compressed iron up from a few GPa to about 140 GPa 
and provided values for a in the range 1.0-1.5 xlO^ri^^ ni^^ (see Ref. [3] and refs. therein; also Ref. [3]). These 
results were in line with earlier estimates of tr ~ 1 x 10^57"^ m"-'^ based on geophysical arguments [5]. A few decades 
later, Secco and Schloessin measured the electrical conductivity of pure solid and liquid Fe at pressures up to few 
GPa and then extrapolated their results to high pressures and temperatures typical of the outer core. Noting the 
effects that the extrapolation to higher parameter values had on the function of the density of states, they derived 
values for a in the range 0.66-0.83 xlO^fi"^ m^^, which agrees well with the results previously obtained by Keeler 
and Mitchell [3] , but is slightly smaller than the value measured by Keeler and Royce [3] . More recently, the electrical 
conductivity of pure iron has been measured in shock compression experiments at high pressures and temperatures 
(up to 208 GPa and 5220 K) by Bi et al. [7]. Its value ranges between 0.76 and 1.45 xlO^fi"^ m"^ by reducing the 
pressure from about 200 to 100 GPa. These findings are in agreement with the result of about 1.0 x 10^17"^ m~^ 
measured by Gomi et al. [5] in analogous experiments performed with a diamond-anvil cell under high static pressure 
(up to 65 GPa) at room temperature. They are also close to the value a = 0.8 x 10^17^^ m^^ estimated by Stacey and 
Anderson by assuming that, for a pure metal, the conductivity is constant on its melting curve. First-principles 
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density functional theory (DFT) calculations of solid iron at low temperature and under pressure have been presented 
by Sha and Cohen [10]. Using the Bloch-Griineisen formula they extrapolated their results up to the high pressures 
and temperatures studied by Bi ct al., finding a slightly larger value in the range 1.1-1.8 x 10^51^^ m^^. Recently 
we have calculated the electrical conductivities of pure liquid iron at the conditions of the Earth's outer core using 
DFT with the Kubo-Greenwood formulation (DFT-KG), obtaining values in the range 1.4-1.6 xlO^rJ"-^ m~^. [TT] An 
independent calculation was also performed by de Koker et al. [12] showing similar values. These values agree well 
with the experimental findings of the older shock wave measurements of Keeler and Royce 0], but they are larger 
than the experimental values reported by Bi et al. above 120 GPa [7^. 

The conductivities of the mixtures with typical core compositions were estimated by Stacey and Anderson [S] to 
be about 0.5 x 10^U~^ m~^ in the case of FeNiSi liquid alloys, by extrapolating measured resistivity values of FeSi 
alloys from Matassov [13] experimental data. Our preliminary results of FeSiO liquid alloys at Earth's core conditions 
indicated only a reduction of about 30% of the conductivity of the mixture compared to that of pure iron [11) . 
suggesting therefore values above 1 x lO^ft^^ m^^. 

Estimates for the thermal conductivity k for liquid iron and liquid iron mixtures ranged between 25 and 60 m^^ 
K~^. [^ irH - fTOj Our recent DFT-KG calculations for pure iron at core conditions provided values significantly higher, 
as expected from the relation between the thermal and the electrical conductivity as encoded in the Wiedemann-Franz 
law [17], which we found to be closely followed [TT] . 

The DFT-KG method has been applied to a wide range of problems, including low pressure systems (C, Na, Al, 
Ar, Ga, Pb) [T5H^ . high pressure (H, H-He mixtures and Al) 24 27J and ultrahigh pressure (water) |28j. Here 
we report detailed DFT-KG calculations for a and k of liquid iron mixtures, Feo.82Sio.10Oo.08 and Feo.79Sio.08Oo. 13, 
representative of Earth's core composition. We find a = 1.1 — 1.3 x 10^fl~^ and k = 100 — 160 W K~^, 
the two extremes in the ranges corresponding to the top and the bottom of the core, respectively. These values are 
in close agreement with extrapolations obtained from recent diamond-anvil-cell experimental measurements of the 
electrical resistivity by Hirose et al. ^9] and by Gomi et al. [SUl 131] (who estimated the thermal conductivity from 
the resistivity using the Wiedemann-Franz law). We find the Wiedemann-Franz law to be closely followed also for 
the iron mixtures, although with a lower value for the Lorenz parameter. 

These new values for the conductivities of the liquid outer core have profound implications for our understanding 
of the Earth's thermal history and the generation of the Earth's magnetic field. 

For the typical Earth's core mixtures investigated here we also show results for the partial radial distribution 
functions, mean square displacements, viscosities and the speeds of sound. 

In Sec. Illjwe describe the techniques used in the calculations. The following section contains our results, starting 



with the pressure-temperature profile in the core in Sec. Ill A the structural properties of the mixtures in Sec. IIIB 



ionic transport properties in Sec jlllCj and electronic transport properties in Sec. |III P] Sec. |IV| includes a discussion 
of the implications of our results for the Earth, and Sec. [Vj contains the conclusions. 



II. TECHNIQUES 

First principles simulations were performed using the VASP code [32 , with the projector augmented wave (PAW) 
method '33', '33 and the Perdew-Wang [3S] functional (PW91). The PAW potential for oxygen, silicon and iron have 
the 25^2^*, 3s^3p^ and Sp^As^icf valence electronic configurations respectively, and the core radii were 0.79, 0.8 and 
1.16 A. To calculate the electrical conductivity we also tested two additional iron PAW potentials, with 4s^3(f and 
3s^sp^4:S^3d7 valence configurations and 1.16 and 0.85 A core radii respectively. The potentials with more semi- 
core states included in valence only give contributions to the optical conductivity at high frequencies, but provide 
the same dc conductivity. For this reason, conductivities have been calculated with the PAW potential having the 
4:S^3(f valence electronic configuration. Single particle orbitals were expanded in plane- waves with a cutoff of 400 eV. 
Electronic levels were occupied according to Fermi-Dirac statistics, with an electronic temperature corresponding to 
the temperature of the system. An efficient extrapolation of the charge density was used to speed up the ab initio 
molecular dynamics simulations [36], which were performed by sampling the Brillouin zone (BZ) with the F point only. 
The temperature was controlled with a Nose thermostat [37] and the time step was set to 1 fs. We ran simulations 
for typically 22 ps, from which we discarded the first ps to allow for equilibration. 

In the adiabatic approximation the electrical conductivity of a liquid can be computed by generating a set of ionic 
configurations {Ri} sampling the relevant pressure-temperature conditions, calculating the conductivity on each of 
these configurations and taking the average. For the conductivities we used the first 5 ps of the simulations to extract 
A'' = 40 configurations {Rj; / = 1, N} equally spaced in time. These N configurations were then used to compute the 
electrical conductivity via the Kubo-Greenwood formula as implemented in VASP by Desjarlais [38]. We performed 
simulations at several thermodynamic states spanning the conditions of the Earth's core, details for the simulations 
used to compute the conductivities will be given in Sec. |IIID| 
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The KG formula for the frequency (w) dependent optical conductivity for a system with Nj ions at positions i?/ 
reads: 

"^^^-■■^'^ = ^^Ml - f^(e,;k)]|(*,,k|V„|vI/,,k)P<5(6,-k - e.M ~ ^) (1) 

i J — 1 Q — 1 

where e and m are the electron charge and mass respectively, h is the Plank's constant divided by 27r, is the volume 
of the simulation cell and n the number of Kohn-Sham states. The a sum runs over the three spatial directions, which 
in a liquid are all equivalent. ^',;,k is the Kohn-Sham wavefunction corresponding to eigenvalue e^^k, and F(ei,k) is the 
Fermi weight. The b function is represented by a Gaussian, with a width chosen to be roughly equal to the average 
spacing between the eigenvalues (0.01 eV for a 157-atom system) weighted by the corresponding change in the Fermi 
fimction [38) . Integration over the Brillouin Zone (BZ) is performed using standard methods 39J, and the frequency 
dependent conductivity is obtained by taking the thermal average: 

aH = (^ak(^;i?7mk)) (2) 

k 

where Vl^(k) is the weighting factor for the point k. The dc conductivity ctq is given by the value of cr(a;) in the limit 
cj — 5- 0. The optical conductivity must obey the sum rule: 

where iVe is the number of electrons in the simulation cell. As reported earlier [40j . to converge the sum rule for 
iron one needs to include states up to ~ 150 eV above the Fermi energy, which means using over 10,000 Kohn-Sham 
states for a typical 157 atoms simulation cell with the 8 valence electron PAW iron potential. However, only states 
near the Fermi energy contribute to the zero frequency limit of the optical conductivity. This is illustrated in Fig. [T] 
where we show as function of w computed using both 10,000 and 1,200 Kohn-Sham states for one configuration 
of liquid iron extracted from the ensemble at p = 328 GPa and T ~ 6350 K. It is clear that both calculations give the 
same dc conductivity, so we decided to spot check the sum rule only in a limited number of cases, and then use 1,200 
Kohn-Sham states for the majority of the calculations. 

In a free electron liquid the electronic part of the thermal conductivity Kq and the electrical conductivity ctq are 
related by the Wiedemann-Franz law, L — kq/uqT, where L is the Lorenz number. In a real liquid the validity of 
the Wiedemann-Franz law is not necessarily expected, and in fact a number of exceptions for metals at near ambient 
conditions are known (see e.g. Kittel [ll]). Here we have directly calculated kq using the Chester-Thellung [42] (CS) 
formulation of the Kubo-Greenwood formula, which reads: 

and kq is the value of k{uj) in the limit cj — > 0. The kinetic coefficients L/„i(w) are given by [13]: 

i^(e,,k)]|(*,,k|V„|vI/,,k)px (5) 
- (5(ej-,k - Ej.k - fi^), 

where /j, is the chemical potential. The implementation of the CS formula in VASP is also due to Desjarlais [381 . 

We checked convergence of the conductivities with respect to the size of the system by performing calculations with 
cubic simulation cells including 67, 157 and 288 atoms and using up to 6 k-points to sample the BZ. We found that 
even the smallest cell sampled with the single k-point (1/4,1/4,1/4) gives results converged to better than 1%. We 
then decided to use 157-atom simulation cells. 



III. RESULTS 
A. Pressure-temperature profile 

The Earth's outer core is convecting, and therefore it is assumed to be in a well mixed state with an adiabatic 
pressure-temperature (p, T) profile. This can be determined by fixing the temperature T at the ICB pressure p = 
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329 GPa and following the line of constant entropy up to the CMB pressure p = 136 GPa. In Fig. [2] we show the p, T 
profiles that we used in the present work, obtained by fixing three possible ICB temperatures: the melting temperature 
of pure iron T = 6350 K [H US] (FERRO), the melting temperature of the mixture Feo.82Sio.10Oo.08 T = 5700 K 
(CORE5700) and the melting temperature of the mixture Feo.79Sio.08Oo.13 (CORE5500) [5^. The mixtures are 
two possible estimates for the composition of the Earth's outer core, which match the Preliminary Reference Earth 
Model (PREM) ICB density jump of 4.5% [47], and the more recent ICB density jump of 6.3% proposed by Masters 
and Gubbins '48', respectively. The adiabats are shown as bands, with the low pressure (solid) edge corresponding to 
the actual DFT-PW91 values, and the high pressure edge to pressures increased by 10 GPa, which is the approximate 
amount by which DFT-PW91 underestimates the pressure of solid iron at Earth's core conditions [44 . 



B. Structure 



In Fig. [3] we show the densities p of the three possible cases mentioned in Sec. |III A on the respective adiabats 



They are shown as bands also in this case, with the same meaning as in Fig. [2] The raw DFT-PW91 density of 
the Feo.82Sio.10Oo. 08 (approximated by a simulation cell containing 129 iron atoms, 16 silicon atoms and 12 oxygen 
atoms) mixture matches the PREM density of the liquid side of the ICB by construction 46 , while the mixture 
Feo.79Sio.08Oo. 13 (125 irons, 12 silicons and 20 oxygens) has a slightly lower ICB density but appears to match quite 
well the CMB density. 

In Fig. |4]we plot the bulk sound velocities as function of pressure for pure iron and for the two mixtures, compared 
with PREM. These are defined as Vb = \l {Ks/ p), where Ks ~ ^V{dp/dV)s is the isentropic bulk modulus, with 
V and S the volume and the entropy of the system, respectively. To compute Ks we fitted the pressures computed 
along the adiabats to a Murnaghan equation of state j49) : 



V 



1 



(6) 



where Kq is the zero pressure bulk modulus and K'^ = [dK / dp)p=Q. Interestingly, the bulk sound velocities of pure 
DFT-PW91 iron are very close to PREM, although the agreement is worsened when the pressure correction of 10 
GPa is applied to the DFT-PW91 calculations. By contrast, the pressure corrected values for the Feo.82Sio.10Oo.08 
mixture sit very close to PREM. As expected, the combination of lower temperatures and densities has the effect of 
increasing the bulk sound velocities. 

In Fig. [5]we show the partial radial distribution functions (rdf) gpeFeif), gsiSi{'f)^ goo{^), 9FeSi{f)^ gFeo{f) and 
9Sio{f)- They are defined so that, by sitting on an atom of the species a, the probability of finding an atom of the 
species fi in the spherical shell (r, r + dr) is iTrr'^npgai3(r)dr, where np is the number density of the species /?. Fig. [5] 
(a) shows the partial rdfs for the mixture Feo.82Sio.10Oo.08 at {p,T) — (329 GPa, 5700 K). In agreement with our 
previous simulations on FeO mixtures [SD] , the present partial rdfs show that the distance between neighbouring iron 
and oxygen atoms (obtained from the position of the first peak of gpeoif) at ~ 1.7 A) is significantly shorter than 
the iron-iron distance (^ 2.1 A) or the oxygen-oxygen distance (~ 2.1 A) . This indicates that oxygen atoms have 
two effective radii, one for the interaction with themselves, and a different one for the interaction with an iron atom. 
The relatively large oxygen-oxygen distance suggests than oxygen atoms effectively repel each other at the typical 
single and double bond oxygen distances (1.47 A and 1.21 A, respectively) [50]. The situation is rather different 
for the silicon-silicon and the iron-silicon distances, as the first peaks of gsiSi{i~) and greSiir) are roughly in the 
same place, showing that iron and silicon atoms have one single effective radius when interacting with each other or 
with themselves, and also that this effective radius is similar for the two atoms. For the silicon-oxygen interaction, 
the position of the first peak of gsio{i~) at ~ 1.6 A is at slightly shorter distance than that of gpeoi^): indicating 
that the silicon-oxygen bond is shorter and probably stronger when compared to the iron-oxygen bond. Although it 
could be suggested that oxygen and silicon in liquid iron may precipitate out as Si02, the simulations provided no 
evidence of any phase separation or departure from a well mixed liquid. In Fig. [s] (b) we show the partial rdfs for 
the same mixture at (p, T) — (134 GPa, 4260 K) (close to CMB conditions), and in Fig. [5] (c,d) those for the mixture 
Feo.79Sio.08Oo.13 at {p,T) = (328 GPa, 5500 K) and (p,T) = (134 GPa, 4112 K), respectively. There are no significant 
differences between the corresponding partial rdfs for the two mixtures, and the only effect of moving from CMB to 
ICB conditions is that of decreasing the height of the .goo(^) peak and increasing the height of the gsioiT) peak, 
which implies that the oxygen-silicon coordination number is slightly increased with pressure. No noticeable pressure 
effect is observed for the iron-oxygen interactions. 



5 



C. Ionic transport 



In this section we describe the ionic dynamical properties of the system, including atomic self-difFusion coefficients 
and viscosities. The self-diffusion coefficient D can be obtained from the asymptotic slope of the time dependent 
mean square displacement (MSD) M{t) in the long time limit t — )■ oo: 



N, 



1 



i=l 



6 t-i-oo dt 



(7) 



where (•) has the meaning of thermal average, which is computed as time average over different origins to along the 
molecular dynamics simulation, ri(t) is the vector position at time t of the i-th atom, and Nj is the number of ions. 
In the inset of Fig. |6] we show the value of M(t) as function of time for the iron atoms in the Feo.82Sio.10Oo.08 mixture 
at ICB conditions. It is clear that after a transient of ~ 0.2 ps the linear behaviour of M{t) is well established. The 
initial part of the transient is due to atoms moving freely before collisions begin to occur, and for this reason the MSD 
increases as the square of time. 

An alternative route to the self-diffusion coefficient is through the Green-Kubo (GK) relations, which relate transport 
coefficients and correlation functions |ST] . The self-difFusion coefficient D is given by the integral of the velocity- velocity 
autocorrelation function (VACF) A{t): 

A{t) = {v,{t + to)-Mto)) 
DA{t) = \j drAir) 
D = lim DA{t) 

(8) 

The function A{t)/A{0) is also plotted in the inset of Fig. |6j where it can be seen that correlations quickly decay to 
zero as soon as the atoms begin to collide with each other. In the main part of Fig. |6] we show the value of M(t) 
divided by 6t, and that of DA{t) as function of time, together with their statistical errors computed by analysing the 
scattering of the MSD and VACF of each atom. The two methods provide the same value for D in the limit of large 
t, with DA[t) converging faster. In Fig. [7]and [Sjwe show the values of the self diffusion coefficients of iron, silicon 
and oxygen for the two mixtures. The values are almost constant along the adiabats, with diffusion slowing down 
only marginally with pressure. Interestingly, iron and silicon have very similar diffusion coefficients, despite a factor 
of two difference in their masses, while oxygen has a diffusion coefficient which is more than two times bigger. No 
significant differences are observed for the two different mixtures. 

To compute the viscosities of the mixtures we used the Green-Kubo relation, which relates the shear viscosity rj 
to the integral of the autocorrelation function of the off-diagonal components of the stress tensor (SACF) Capit) = 
{Pa/sit + to) Pai3{to)), where Pais is an off-diagonal component of the stress tensor, with a and (3 indicating Cartesian 
components. There are five independent components of the traceless stress tensor: P^y, Py^, Pzx, ^{Pxx — Pyy) and 
\{Pyy ~ Pzz) In a liquid these five components are equivalent, and they can all be used to improve the statistical 
accuracy of the viscosity integral. The shear viscosity rj is then obtained from: 



V 



kbT Jo 

rj — lim r/(i) 



(9) 

where fcs is the Boltzmann constant and S{t) = l[Cxy{t) + Cyz{t) + Czx + ^{Cxx{t) ~ Cyy{t)) + \{Cyy{t) — Czz{t))\- 
In Fig. [9] we plot S{t)/S{Q) and rj{t) for the Feo. ,82Sio.ioOo.o8 mixture at ICB conditions. To obtain the shear viscosity 
T) we need to integrate S{t) to i = 00, however, it is clear that the S{t) has decayed to zero after ^ 0.2 ps, after 
which one only integrates statistical noise. For this reason, we decided to stop the integration dX t — 0.2 ps. The 
viscosity integral rjit) is also plotted in Fig. [oj together with its error bar estimated by the scattering of the five 
independent components of the traceless SACF. To improve the statistics on the estimate of the statistical error on 



7y(t) we split the simulations in two independent chunks. The viscosities are plotted in Fig. 10 for the two mixtures on 
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the corresponding CORE5700 and CORE5500 adiabats. Their values range between ~ 7 mPa s at the CMB and ~ 12 
mPa s at the ICB, with httle difference between the viscosities of the two mixtures. We also report the viscosities of 
pure iron on the FERRO adiabat, which are consistent with those reported before [551 [51] . and they have roughly the 
same values as those of the mixtures. 

To conclude this section we also report the ionic component of the thermal conductivity of pure iron, Kion, computed 
using a simple pair potential, which describes the energetics and the structural and dynamical properties of iron very 
accurately [54]. The pair potential has the form U{ri, . . . , rjvj) = 4:J2f=i i<i{'^/\^i~''^j\)" ^ where {r^} are the Cartesian 



coordinates of the atoms, a = 5.86, and a = 1.77 A. The Kubo-Green formula for the ionic thermal conductivity 
reads: 

1 



OO 



mkeT^ JO 

where the microscopic heat current is given by: 

Ni Ni 



dmt + to) -iito)) (10) 



jW = Ev,e, + i r.,(F,, -v,) (11) 

with Wi the velocity of atom i, Vij the vector distance between atom i and atom j, and F^j the force on atom i due 
to atom j from the pair potential. The on-site energy is: 



1 . Ni 



with m being the mass of the atoms. The calculated values range between 2.5 and 4 W m -'^ K which compared 
to the electronic contribution to the thermal conductivity (see next section) are completely negligible. 



D. Electronic transport 

The electrical and thermal conductivity of pure iron on the three adiabats shown in Fig. [TT] have been reported 
before 11^, where we also mentioned that preliminary results showed that light impurities have the effect of reducing 



the conductivities of pure iron by 30%. Here in Fig. 12 we report the full calculations on the CORE5500 and the 



CORE5700 adiabats for the mixtures. It is clear that indeed the effect of the light impurities is that of reducing 
the conductivities of pure iron by the reported amount of 30%, and that the slightly different composition of the 
two mixtures has roughly the same reduction effect. The values for the electrical and thermal conductivities are in 
the range 1.1-1.3 xlO^il"^ m~^ and 100-160 W m^^ K~^, respectively, with the low/high values corresponding to 
CMB/ICB pressure-temperature conditions. 

The Lorcnz parameter is roughly constant on all three adiabats, indicating that the Wiedemann-Franz law is valid 
throughout the core. For pure iron, the Lorenz parameter varies between 2.47 xlO^* W ft K^^ and 2.51 xl0^8 W 

K~^, only slightly higher than its ideal value of 2.44 xlO^^ W H. K^^, while for the mixtures it is reduced in the 
range 2.17-2.24 xlO"** W n K'^. 

The values we find for the conductivities and the Lorenz number are in broad agreement with those recently reported 
by de Koker et al. [H] . Our electrical conductivities are also in fairly good agreement with the experimental data for 
FeSi up to 140 GPa of Matassov [T3j and for Feo.940 up to 155 GPa of Knittle et al. [SS], who reported values in the 
range 1.0-1.2 xl0^fl~^ m~^. Our thermal conductivities are also in agreement with recent experimental findings of 
Hirose et al. [5S] who reported values for the top of the outer core in the range 90-130 W m~^ K~^. 

All results are summarised in Table I, II and III. 



IV. IMPLICATIONS FOR THE EARTH 



The estimates of thermal and electrical conductivity in Tables [TT] and [TTT] are 2 — 3 times higher than those currently 
used in the geophysical literature [e.g. 156 1 [57 ] . The high thermal conductivity in particular has significant implications 
for the evolution of the core and the dynamo process generating the Earth's magnetic field. The convective motions 
in the outer core that are responsible for the Earth's dynamo are driven by a combination of thermal and chemical 
buoyancy sources. The strength of thermal driving is measured by the amount of excess heat that cannot be conducted 
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down the adiabatic gradient; higher thermal conductivity increases adiabatic conduction and therefore decreases the 
effectiveness of thermal buoyancy relative to chemical buoyancy. Maintaining the same magnetic field with less 
available thermal buoyancy requires a faster core cooling rate or a higher concentration of radiogenic elements in the 
core or a combination of the two. Moreover, a faster cooling rate implies that the inner core, which is already thought 
to be a relatively young feature of the Earth (age ~ IGa, |SH]), is even younger. 

The reduction in thermal buoyancy due to increased adiabatic conduction is so great that parts of the liquid core 
are very likely to be subcritical to thermal convection. Chemical convection, driven by release of light elements at the 
boundary of the solid inner core on freezing, may still be able to stir these regions, heat being transported downwards 
to maintain the adiabatic gradient [SH] . Near the top of the core chemical convection weakens because of the barrier of 
the core-mantle boundary (CMB); here the liquid is likely to be density stratified with little or no vertical motion at 
all [TT] . Such a stable layer could be detected observationally, either by seismology or by its effect on the geomagnetic 
field. 

It has been suggested that the heat convected away from the inner core may vary so much from place to place 
that the surface may be melting in some places, freezing in others, with a net growth from freezing over the whole 
surface [60 . The increase in thermal conductivity with depth (Tables [TT| and III I increases the heat conducted away 
from the ICB, making the chance of melting rather less likely. It also reduces the vigour of thermal convection deep 
in the core. 

The main consequence of a higher electrical conductivity is to lengthen the ohmic diffusion time, the main time 
scale used when interpreting geomagnetic and paleomagnetic phenomena. Doubling or trebling the time scale affects 
virtually all interpretations to a certain extent: the dipole decay time increases from 25 kyr to 75 kyr. For example, it 
improves the " frozen-flux" assumption, commonly used to interpret recent geomagnetic secular variation j61j , changes 
on time scales of decades-to-centuries; a polarity reversal that can take anything from 1 to 10 kyr now appears fast 
on the diffusion time scale. The magnetic Reynolds number for the dynamo is also larger, making dynamo action 
possible with lower flow speeds. Lastly, the inner core has been thought to provide a stabilising influence because its 
ohmic diffusion time (5,000 years) is longer than the advective time in the outer core (500 years) [62] ■ The higher 
electrical conductivity reported in this paper increases the diffusion time considerably, making the stabilising effect 
even stronger. 

The ab-initio results also have implications for numerical models of the geodynamo. The equations governing the 
geodynamo process are usually cast into nondimensional form; it is well known that fundamental aspects of the 
solutions depend critically on the values of the dimensionless parameters. Common dimensionless parameters are 
the thermal Prandtl number Prx — v/Dt, the ratio of viscous and thermal diffusion, the magnetic Prandtl number 
Pm = v/Db, the ratio of viscous and magnetic diffusion, and the chemical Prandtl number Prx = v/D-y^, the ratio of 
viscous diffusion and mass diffusion for species X =0,Si. Here v = r}/ p is the kinematic viscosity, Dt — hq/ pCp is the 
thermal diffusivity where Cp is the specific heat at constant pressure, and Db = 1/(mo'7o) is the magnetic diffusivity 
where /lo is the permeability of free space. 

In Table IV we list values for the Prandtl numbers using the results of our ab initio calculations (Tables |ll] and 



nil at the CMB (P=134 GPa) and ICB together with errors. Pm is slightly higher than a previous estimate of 
6 X 10~^ [55], but this is still too low to be achievable in a simulation in the near future. Values of Ptt are lower than 
the value Ptt = 0.1 commonly assumed for liquid metals and can certainly be accommodated in current geodynamo 
models. The values of Prgi and Pro differ significantly and are larger than previous estimates [63| . The Lewis number 
Le — Pvq/ Ptt, an important parameter in dynamo models driven by both thermal and compositional buoyancy, is 
found to be 0(10^). 



V. CONCLUSIONS 



In this work we have studied the transport properties for two liquid silicon-oxygen-iron mixtures, i.e., 
Feo.79Sio.08Oo.13 and Feo.82Sio.10Oo.08> at Earth's core conditions using DFT-KG ab-initio theoretical calculations. 
We find that both the thermal and electrical conductivity of the mixtures are higher than previous estimates, the 
former being in the range of 100-160 W K-i and the latter in the range of 1.1-1.3 xlO^fJ ^ m ^. The validity 
of the Wiedemann-Franz law is found to be satisfied quite accurately for both iron mixtures at core conditions, the 
Lorenz parameter being roughly constant within the range 2.17-2.24 xl0^8 W K^^, only slightly higher than the 
ideal value of 2.44 xl0"8 W ft K~^. The values we find for the conductivities are 2 to 3 times higher than those 
which have been used until now and have profound implications for our understanding of the outer core evolution and 
the geodynamo. In particular, the inner core is now younger than previously thought, the top of the core is very likely 
to be thermally stably stratified, while the possibility that the inner core is partially melting is less likely. Finally, 
we list values for the dimensionless input parameters used in geodynamo simulations, calculated directly from the ab 
initio calculations. 
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10.48 
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10.36 


0.50(1) 
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10.52 


7.93 


0.59(1) 


6.2(4) 


1.338(5) 


154(1) 


2.49(2) 



TABLE I: Pressure (P), temperature (T), density (p), bulk sound velocities («;,), iron diffusion coefficient (Dpe), viscosity {rj), 
electrical conductivity {(To), thermal conductivity (/to) and Lorenz parameter (L) for pure liquid iron on the FERRO adiabat. 
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TABLE IL Pressure (P), temperature (T), density (p), bulk sound velocities (vb), diffusion coefficients for iron (Dfc), silicon 
(Dsi) and oxygen (Do), viscosity (r;), electrical conductivity (ao), thermal conductivity (/to) and Lorenz parameter (L) for the 
Feo.82Sio.10Oo.08 liquid mixture on the CORE5700 adiabat. 
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TABLE III: Pressure (P), temperature (T), density (p), bulk sound velocities (vb), diffusion coefficients for iron (Dpe), silicon 
(Dsi) and oxygen (Do), viscosity (r?), electrical conductivity (cro), thermal conductivity (/to) and Lorenz parameter (L) for the 
Feo.79Sio.08Oo. 13 liquid mixture on the CORE5500 adiabat. 
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328 0.063 (0.004) 
134 0.048 (0.005) 


1.7 (0.12) 
0.9 (0.11) 


286 (27) 118 (12) 
147 (24) 52 (9) 
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1.5 (0.13) 
1.0 (0.09) 
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98 (14) 
52 (7) 


1.8 (0.03) 
1.5 (0.02) 



TABLE IV: Values of the thermal Prandtl number (Prr), magnetic Prandtl number (Pm), and chemical Prandtl numbers for 
Silicon (Prsi) and Oxygen (Pro). Also reported is the thermal diffusivity Dt- Top section is for the inner-core density jump 
of Masters and Gubbins 48 ; bottom section is for the PREM density jump [l^. Error estimates are in brackets. 
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FIG. 1: Optical conductivity u{(jj) of liquid iron as function of energy computed using Nks =10,000 and Nks =1,200 
Kohn-Sham states for one configuration extracted from the ensemble at p = 328 GPa and T = 6350 K. 
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FIG. 2: Adiabatic temperature profiles as function of pressure anchored by three possible ICB temperatures Ticb- Black curve 
(FERRO) corresponds to Ticb = 6350 K. Red curve (CORE5700) to Ticb = 5700 K, and green curve to Ticb = 5500 K. Low 
pressure edges of the bands (solid line) correspond to the raw DFT-PW91 pressures, high pressure edges to pressure corrected 
by 10 GPa, which is the typical DFT-PW91 underestimate of the pressure for iron at core conditions (see text for details). 
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FIG. 3: Density profiles for pure iron (black), the Feo.82Sio.10Oo.08 mixture (red), and the Feo.79Sio.08Oo. 13 mixture (green) 
on the respective FERRO, CORE5700 and CORE5500 pressure-temperature profiles displayed in Fig. [2] Bands have the same 
meaning as in Fig. [2] PREM density profile [47 is shown with blue line. 
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FIG. 4: Bulk sound velocity for pure iron (black), the Feo.82Sio.10Oo.08 mixture (red), and the Feo.79Sio.08Oo. 13 mixture (green) 
on the respective FERRO, CORE5700 and CORE5500 pressure-temperature profiles displayed in Fig. [2] Bands have the same 
meaning as m Fig. [2] The PREM bulk sound velocities [47] are shown in blue. 
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FIG. 6: Self-diffusion coefficient as a function of time obtained from tlie time dependent mean square displacement M{t) (red) 
and the velocity- velocity autocorrelation function Da (i) (black) (see text details). Also shown in the inset the value of Al{t) 
and A{t)/A{0) as a function of time. Results correspond to the iron atoms in the Feo.82Sio.10Oo.08 mixture at ICB conditions. 
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FIG. 7: Self diffusion coefficients of Fe, Si and O atoms vs pressure for the Feo.82Sio.10Oo.08 mixture. Temperatures are 
determined from tfie CORE5700 adiabat at the corresponding pressures. 
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FIG. 8: Self diffusion coefficients of Fe, Si and O atoms vs pressure for the Feo.79Sio.08Oo. 13 mixture. Temperatures are 
determined from the CORE55Q0 adiabat at the corresponding pressures. 
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FIG. 9: Average over the five independent components of tlie autocorrelation function of tlie traceless stress tensor S{t), 
normalised by dividing by ^(O). Also shown is the viscosity integral ri{t) with its statistical error. Results correspond to the 
Feo.s2Sio.10Oo. 08 mixture at ICB conditions (see text for details). 
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FIG. 10: Viscosity as a function of pressure for the pure iron on the FERRO adiabat, and for the mixtures Feo.82Sio.10Oo.O8 
and Feo.79Sio.08Oo. 13 on the CORE5700 and CORE5500 adiabats, respectively. 
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FIG. 11: Electrical (a) and thermal (b) conductivity of liquid iron at Earth's core conditions, computed on the FERRO (black), 
CORE5700 (red) and CORE5500 (blue) adiabats. Lines are quadratic fits to the first principle raw data (symbols). Error bars 
(2 s.d.) are estimated from the scattering of the data obtained from 40 statistical independent configurations. Results are 
obtained with cells including 157 atoms and the single k-point (1/4,1/4,1/4), which are sufficient to obtain convergence within 
less than 1%. 
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FIG. 12: Electrical (a) and thermal (b) conductivity of liquid Feo.82Sio.10Oo.08 (red) and liquid Feo.79Sio.08Oo. 13 (green) 
mixtures on the CORE5700 and CORE5500 adiabats, respectively. Error bars (2 s.d.) are estimated from the scattering of 
the data obtained from 40 statistical independent configurations. Results are obtained with cells including 157 atoms and the 
single k-point (1/4,1/4,1/4), which are sufficient to obtain convergence within less than 1%. 



